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ABSTRACT 

We discuss the propagation of electromagnetic waves on a rectangular lattice of po- 
larizable point dipoles. For wavelengths long compared to the lattice spacing, we obtain 
the dispersion relation in terms of the lattice spacing and the dipole polarizabilities. We 
also obtain the dipole polarizabilities required for the lattice to have the same disper- 
sion relation as a continuum medium of given refractive index m; our result differs from 
previous work by Draine & Goodman (1993). Our new prescription can be used to 
assign dipole polarizabilities when the discrete dipole approximation is used to study 
scattering by finite targets. Results are shown for selected cases. 

1. Introduction 

The discrete dipole approximation (DDA) provides a flexible and general method to calculate 
scattering and absorption of light by objects of arbitrary geometry (Draine 1988; Draine & Flatau 
1994). The approximation consists of replacing the continuum target of interest by an array of 
polarizable points, which acquire oscillating electric dipole moments in response to the electric field 
due to the incident wave plus all of the other dipoles. Broadly speaking, there are two criteria 
determining the accuracy of the approximation: 



• the interdipole separation d should be small compared to the wavelength of the radiation 
in the material {\m\kd < 1, where m is the complex refractive index, and k = to/c is the 
wavenumber in vacuo); 

• the interdipole separation d should be small enough to resolve structural dimensions in the 
target. 

With modern workstations, it is now feasible to carry out DDA calculations on targets containing 
up to ~ 10 6 dipoles (Draine &: Flatau 1994). In addition to calculation of scattering and absorp- 
tion cross sections, the DDA has recently been applied to computation of forces and torques on 
illuminated particles (Draine & Weingartner 1996, 1997) 
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If the dipoles are located on a cubic lattice, then in the limit where the interdipole separation 
kd — > 0, the familiar Clausius-Mossotti relation (see, e.g, Jackson (1962)) can be used to determine 
the choice of dipole polarizabilities required so that the dipole array will approximate a continuum 
target with dielectric constant e. Draine (1988) showed how this estimate for the dipole polar- 
izabilities should be modified to include radiative reaction corrections, and Draine & Goodman 
(1993) derived the 0[(kd) 2 ] corrections required so that an infinite cubic lattice would have the 
same dispersion relation as a continuum of given dielectric constant. 

Nearly all DDA calculations to date have assumed the dipoles to be located on a cubic lattice. 
If, instead, a rectangular lattice is used, it will still be possible to apply FFT techniques to the 
discrete dipole approximation, in essentially the same way as has been done for a cubic lattice 
(Goodman, Draine & Flatau 1991). However, the ability to use different lattice constants in 
different directions might be useful in representing certain target geometries. 

The objective of the present report is to obtain the dispersion relation for propagation of elec- 
tromagnetic waves on a rectangular lattice of polarizable points. This dispersion relation will be 
expanded in powers of (kd), where d is the characteristic interdipole separation. For a lattice of 
specified dipole polarizabilities, this will allow us to determine the dispersion relation for electro- 
magnetic waves propagating on the lattice when (kd) 2 <C 1. Alternatively, if we require the lattice 
to propagate waves with a particular dispersion relation, inversion of the lattice dispersion rela- 
tion will provide a prescription for assigning dipole polarizabilities when seeking to approximate a 
continuum material with a rectangular lattice of polarizable points. 



2. Mode Equation for a Rectangular Lattice 

The problem of wave propagation on an infinite polarizable cubic lattice of point dipoles 
has been solved previously by Draine & Goodman (1993, herafter DG93). Here we describe the 
generalization of this problem to rectangular lattices. 

Consider an infinite rectangular lattice with lattice sites at 

x n = (n 1 d 1 ,n 2 d 2 ,n 3 d 3 ) , (1) 

where the rij are integers and the lattice constants di are, in general, all different. The density n of 
dipoles is just 

n = d^dS ' (2) 

and it will be convenient to define the characteristic lattice length scale 

d=(d 1 d 2 d 3 ) 1 / 3 . (3) 

Since the lattice is anisotropic, the polar izability a. of the dipoles located at lattice sites is a tensor, 
i.e., the polarization is 

P(x) = a ■ E(x) . (4) 

Thus the polarization vector is not, in general, parallel to the electric field or to the vector potential. 
From charge conservation we have 

V-J = 0, (5) 
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or, since 

we obtain the transversality condition 

V • P = . (7) 
We assume that the dipole moment at the lattice site x n is 

P n (t) = P (0)e ik - Xn ^. (8) 

With the Lorentz gauge condition 

1 <9$ 

VA + -_ = (9) 
and equation (6), the vector potential A(x, t) satisfies the wave equation 

V 2 A + ^A = (10) 

& c 

lr/ ^P n ^(x-x n ) (11) 



c 

n 



^P e^V<5 3 (x-x n ) . (12) 
f — ^ 



Consider now the unit cell centered on x = 0. Let A ot her and Mother be the potentials in this 
region due to all dipoles except the dipole at x = 0. Thus 

B other (*,t) = --^^ ~ Mother = ^A!(x,t) (13) 

where 

c 2 

Ai (x, t) = A other + — V(V • A other ) . (14) 



Thus 



P (0) = — a • Ai(x = 0,t = 0) . (15) 



The vector potential can be written 

A(x, t) = e *-*-*"* ^ a n e iq ' x (16) 



q = q(n)=27r — ,— ,— (17) 



where 

. di cfe d 3 

define a lattice which is reciprocal (e.g., Ashcroft & Mermin (1976)) to the original rectangular 
lattice on which the dipoles are located. Equations (12) and (16), and the identity 

E *"<*-*■>- 55* (18) 
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yield 



1 



C 2 diG?2^3 



a • Ai(0,0) 



2-1 - 1 



|k + q| 2 --o 



The field A se ^(x, t) due to the dipole at x = is 



A se// (x,f) 



-»w e 



iuir/c Jl 



-=-a-Ai(0,0)- 



O iior/c 



-iuit 



With the identity 



we obtain 



*-other 



4tt^ 2 1_ 

c 2 did 2 <^3 



2vr 2 



«k-x— iuit 



S Sk 'w 



ik'-x 



+ k| 2 - to 2 /c 2 



(19) 

(20) 
(21) 



a- Ai(0,0) 



did 2 d 3 



E 



3l/ 



k + q| 2 — lo 2 /c 2 
„ik'x 



8vr 3 y |k + k'| 2 -lo 2 /c 2 
We seek Ai. Substituting (22) into (14), we obtain after some algebra: 



(22) 



3 3 



47T 



i(k+q)-x 



8vr 3 



j(k+k')-x 



|k + k'| 2 -co 2 /c 2 _ 



k + q| 2 — {J 1 /c 2 
^«^i, 3 (0,0). 



(23) 



For x = and t = and given values of otij, equation (23) is the mode equation, i.e., it allows 
determination of the dispersion relation k = k(cj). However, in the present case, we know the 
dispersion relation, and wish to determine values of ay such that the known dispersion relation 
and equation (23) are satisfied simultaneously. Another relation that must be satisfied follows from 
(7), (8), (15), and (18): 

k-Po(O) =k a- Ai(0,0) = 0. (24) 



Hence, if e is a unit vector in the direction of Ai, then 

3 3 

EE*" ai i e i = • 

j=l 1=1 

At this point it is convenient to define the following dimensionless quantities: 

2vr 2tt 

q = 



d (d^dz) 1 / 3 



V = 



qc ' 



= 



(25) 

(26) 
(27) 
(28) 
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4-7T 
did 2 d 3 

q , d d d 
Q = - = ni7> n 27,n37 , (30) 
q di d 2 d 3 

The vectors Q define a reciprocal lattice (e.g., Ashcroft & Mermin (1976)). 

Using these dimensionless quantities, the mode equation (23) may be written 

3 3 

Y Mil TV e i = ei ' ( 31 ) 

3=1 1=1 

where 

13 V |/3 + Q| 2 -z/ 2 J \(3 + 0\ 2 -v 2 ' 1 ' 

and equation (25) becomes 

3 3 

EE^«ej=0- (33) 
i=i j=i 

It should be noted that the mode equation for a cubic lattice has the same form as equation (31). 
The matrix elements Mjj in equation (32) also have the same form as the matrix elements for 
a cubic lattice but with the vector n replaced by Q for a rectangular lattice. Since the matrix 
elements for a cubic lattice have been calculated previously, it is convenient to rewrite (32) as 

M l3 =M% + T lJ , (34) 
where M^j is the matrix element for a cubic lattice, and Tjj is 

rp _ \- v 2 $ij - (Pi + Qi)(Pj + Qj) \- v 2 8jj - (Pi + m)(Pj + rij) 

13 1/3+ Q|2 _ ^2 2^ \p + n \2_„2 ' 1 j 



3. Dispersion Relation in the Long- Wavelength Limit 

For a given polarizability tensor 7, the mode equation (31) allows one to determine the dis- 
persion relation (3 = (3(v). In the present case, the dispersion relation is known: 

(3 = mva. , (36) 

where m = e 1 / 2 is the (complex) refractive index and a is a unit vector in the direction of propaga- 
tion, which is fixed. From equation (33) for a rectangular lattice, the polarization rather than the 
vector potential is perpendicular to the direction of propagation, i.e., 



3 3 

Y Y ai ^ e i = • 

i=l 3=1 



(37) 
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Therefore, the mode equation (31) can be used to determine dipole polarizabilities 7^ such that 
equation (37) is satisfied. 



We seek the dipole polarizabilities in the limit v 2 <C 1. In this limit, one can write 

Mg = M^+M^+0^) , 
T l3 =T^+T^ + 0^) , 
e = e (°) + e« + 0(z/ 4 ) , 



(38) 

(39) 
(40) 
(41) 



where superscripts and 1 designate the zero frequency limit and the leading order corrections, 
respectively. 

The low-frequency expansion of the matrix element Mjj for a cubic lattice was given by DG93 
and is: 



Mr. 



C(0) 



m 2 + 2 



m 



v 3(m 2 -l) 13 m 2 -l 



a; a 



l u J 1 



= u 2 5ij [ci + m 2 c 2 (l - 3a 2 )] + m 2 ^ 2 (l - 5ij)c 3 a,iaj 



4n 2 i 



where 



(42) 
(43) 



d = -5.9424219... 
c 2 = 0.5178819... , 
c 3 = 4.0069747... . 



From equation (35), one finds 
where 



r (o) 



5ijR (i) , 



IQI 



(44) 



(45) 



Note that Yli=i Ro(i) = 0, so that there are only two independent sums which require evaluation. 
Numerical evaluation is discussed in the Appendix. 



ij 



{1) - m 2 u 2 a iaj [R x - 2R 2 (i) - 2R 2 (j) + 8R 3 (i,j)] 



u 2 5„ 



Kl 



Ri + (m 2 - l)R 2 (i) + 8m 2 a 2 R 3 (i,i) - Am 2 ^ a 2 R 3 (i, I) 



1=1 



where 



(46) 



(47) 
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IQI 

n 2 n 2 2 2 
n i i n 1^*1 

In the limit of cubic lattice, Qi — > nj, and all the sums -R n in equations (45) - (49) vanish. For 
a general rectangular lattice, these sums must be evaluated for each set of new lattice constants. 
Care must be used in evaluation of these sums as described in the next section. We note that R\ 
and i?2 may be obtained from R3: 

3 3 

i=l j=l 
3 

R 2 (i) = J2Mh3) • (51) 
Since Rs(i,j) = Rs(j,i), six independent quantities suffice to determine R\, R2, and R3. 



3.1. Static Limit 

In the static limit the polarizability tensor 7^ is diagonal, and the mode equation and transver- 
sality condition become: 



EK^+^H^ef =ef\ (52) 



^^^ = 0. (53) 



=1 



Substituting (42) and (44) into (52) and taking into account condition (53) one gets 



m 2 + 2 . 

+ Rq{i) 



3(m 2 - 1) 



(0) (0) (0) 
7ii e ; = e, v • (54) 



For e- -* 7^ 0, the static dipole polarizability for a rectangular lattice is 



(0) _ 1CM , , 

lu ~l+7CMRo(i)' (55) 
where the Clausius-Mossotti polarizability is 

_ 3(m 2 - 1) 

1CM = 2 i o • 56 ) 

It is evident from equation (55) that the static polarizability for a cubic lattice reduces to the 
Clausius-Mossotti result. 

If ej°) = in one of the directions, then the polarizability in this direction is undetermined by 

equation (54). However, since the result (55) is valid for values of arbitrarily close to zero, we 
take (55) to be the solution for all directions of the electric field subject to equation (53). 
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3.2. Finite Wavelength Corrections 



We now proceed to obtain the leading corrections for small v. The first order mode equation 
and transversality condition are: 



£ £ 4 0) + £ £ < K" 4 0) + 7,™ 4" 

j=i ;=i i=i ;=i 



= e„ (1) 



(57) 



EX> [■ 

i=i i=i 

From equations (42), (44), (54), and (58) one gets 

3 



7 « e (0) + 7 (0) e (l) 
3 ^ J V 3 



0. 



3 3 

EE^ 



(0) 



Jl) (0) (0) (1) 

7 «i e i e i 



E 



m 2 + 2 



^ 3(m 2 - 1) 
(o) e i • 



7 « e (. ) +7 (°) e w 



3 

JO) j_ Ji) 

. 7 (0) 
j=l h 



Therefore, equation (57) becomes 



3 

E 



(0) ^ y '33 

"in 



ef =0 



It is convenient to write M-^ in the form 



= v 2 5ijLn + m 2 u 2 aiajKij 



where, from (43) and (46) 



La = c\ + m z C2(l — 3a 2 ) — m 2 c 3 a 2 



4ir 2 i 



(58) 



(59) 



(60) 



(61) 



- Ri - im 2 - l)R 2 (i) - 8m 2 a 2 R 3 (i,i) + 4m 2 ^ a 2 i? 3 (M) (62) 

l=i 

Kij = c 3 + R 1 - 2R 2 (i) - 2R 2 {j) + 8R 3 (i,j) . (63) 

We now observe that one can have the square bracket in equation (60) vanish for each j if the 
first-order correction to the polarizability is 



and 



7 w = - 7 (°U°)m (1) 

hj hi Ijj 13 

= ~ lj2 "lii } lfj [Lutiij +m 2 a i a j K ij ] . 



(64) 
(65) 
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Equation (65) makes it explicit that the frequency-dependent term in the dipole polarizability 
depends on the direction of propagation of the wave. 

We observe, however, that this solution is not unique. In fact, the transversality condition (58) 
permits the more general solution 

llf = -v 2 lf ifj [LuSij + m 2 a iaj (K tJ + c 4 )] (66) 

where c 4 could, in principle, be a function of m 2 , the direction of propagation, or the polarization. 
We note, however, that if C4 7^ then eq. (60) is no longer satisfied term by term, but only after 
summation. 

Since we appear to have freedom in the choice of C4, numerical experiments using different 
choices for C4 could be used to see what value of C4 appears to be optimal for scattering calculations 
using the discrete dipole approximation.. 



3.3. Special Case of a Cubic Latttice 



Propagation of electromagnetic waves on a cubic lattice has been discussed by DG93. In their 
analysis it was assumed that if the dielectric material is isotropic, the polarizability 7^ = 7^ +7^ 



no 



can be taken to be isotropic: 7^ = j^5ij and 7^ = ^^5^. With this assumption, contraction 
of eq. (60) with leads to 



(i) 



(o) 
^33 



v Si 



An 1 



c\ + m C2 — m (3c2 + 03)8 — iu 



(67) 



(68) 



j'=i 



as in eq.(4.14) of Draine & Goodman (1993). However, we have seen above that for a cubic lattice 
(for which Kij = C3), 7-^ can be made diagonal but not isotropic: if we choose C4 = — C3 so that 



1 ij 

7^ is diagonal, we obtain 



hj 



(0) 
^33 



4vr 2 



ci + m c-2-m (3c 2 + c 3 )a ? - — iu 



(69) 



which differs from the result of Draine & Goodman (1993) since a 2 - ^ S, except for the special case 
of propagation in the (1,1,1) direction, for which af = S = 1/3. 



4. DDA Scattering Tests 



We now compare discrete dipole approximation (DDA) calculations of scattering and absorp- 
tion using (1) the polarizability prescription (67) from DG93 and (2) the present result (69). We 
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Fig. 1. — DDA calculations of scattering and absorption by a sphere with refractive index m = 1.33 + 0. OH. Upper 
panel shows Q SC a, Qabs, and Q cxt = Q aC a + Qabs calculated using Mie theory (solid curves) and using the DDA with 
N = 17904 and the present polarizabilities. Middle and lower panels show fractional errors for Q sca and Q a bs- Circles: 
using dipole polarizabilities from this paper (eq.(69). Triangles: using dipole polarizabilities from Draine & Goodman 
(1993); 



consider spheres where we can use Mie theory to easily calculate the exact result for comparison. 
Let C sca and C a b s be the cross sections for scattering and absorption, and let Q sca = C scs ,/ira 2 and 
Q a bs = Cabs/ira 2 be the dimensionless scattering and absorption efficiency factors for a sphere of 
radius a. 

Figure 1 shows Q sca and Q a hs calculated using the DDA with a target consisting of N = 17904 
dipoles, with radius a = a e ff = (3A^/4-7r) 1 / 3 (i. The DDA calculations used a cubic lattice, and 
were carried out using the public domain code DDSCAT.6.1 (Draine & Flatau 2004). Calculations 
were done for values of the scattering parameters x = 2-7ra/A satisfying \m\kd < 0.8, as this 
appears to be a good validity criterion for DDA calculations. The upper panel shows scattering and 
absorption efficiencies calculated with Mie theory and using the DDA with the present prescription 
for polarizabilities a. The middle and lower panels show fractional errors [Q (DDA) /Q (Mie) — 1] of 
the DDA scattering and absorption cross sections. 

Because the DDA target (dipoles on a cubic lattice) is not rotationally symmetric, the scat- 
tering problem was solved for 12 different incident directions, in each case for two orthogonal 
polarizations. The cross sections shown in Fig. 1 are averages over the 12 incident directions and 
two incident polarizations. 



-11- 




In the limit x — > 0, identical results are obtained for both polarizability prescriptions; this is 
expected since they differ only in terms of 0(\m\ 2 (kd) 2 ). Note that even at x = (i.e., kd = 0) - 
for which the polarizabilities are exactly given by the Clausius-Mossotti polarizabilities - the DDA 
calculation has a nonzero error even as x — * and kd — * 0. This is because (1) the dipole array only 
approximates a spherical geometry and (2) when the "lattice dispersion relation" polarizabilities 
are used, the DDA does not accurately model the "screening" of the target material close to the 
surface by the polarization at the target surface. For modest dielectric functions, and sufficient 
numbers of dipoles, the errors are small - for the case shown in Figure 1 the fractional errors in 
the x — > limit are < 0.3%. 

As x increases, the two prescriptions result in differences in computed values of Q sca and Q a bs> 
although the differences are small, and comparable in magnitude to the errors in the x — > limit. 
In Figure 1 our new prescription appears to result in small improvements in accuracy (relative to 
DG93) for Q sca for \m\kd < 0.8, and likewise for Q a b s for \m\kd < 0.5. For 0.5 < \m\kd < 0.8, 
however, the DG93 prescription appears to give more accurate results for <5 a bs- Note, however, 
that for both polarizability prescriptions the fractional errors in Q sca and Q a b s for [m|A;d « 0.5 are 
of the same order as the fractional errors for kd — > 0, so it is not clear that the errors should be 
attributed to the lattice dispersion relation - the errors in Q sca and Q a b s may instead be associated 
with errors in the computed polarization field near the surfaces, arising from the fact that in the 
DDA the dipoles in the surface layer are imperfectly "screened" from the external field by their 
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neighbors. 
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Fig. 3. — As in Figure 1, but for m = 2 + i. 

For \m\kd < 0.8, the errors in Q S ca using the new polarizabilities are of order ~ 1% for 
m = 1.33 + 0.01i (Figure 1), m = 1.7 + O.li (Figure 2) and m = 2 + i (Figure 3), and of order ~ 5% 
for m = 3 + 4i (Figure 4). In all four cases, for wavelengths satisfying \m\kd < 0.8, the magnitude 
of the error in Q &CA is of the same order as the error at zero frequency, suggesting that it is mainly 
due to the surface. 

The errors in Q a b s tend to be somewhat larger than for Q sca ,, especially for large values of 
Im(m). For m = 3 + Ai (Figure 4), the new polarizabilities result in a fractional error of ~ 50% for 
x = 2.6 and N = 17904 dipoles. 

Figure 5 shows results calculated for m = 1.7+0.1i spheres represented by N = 140408 dipoles 
- 8 times as many as in Fig. ??. Note that the fractional error for x — > is now only 0.5% - half 
as large as for N = 17904. This can be ascribed to the fact that the fraction of the dipoles in the 
surface monolayer was reduced by a factor of 2. Note also that the errors for finite x are also about 
a factor of 2 smaller than for N = 17904. These results appear to confirm the suspicion that the 
errors arise primarily from the surface layer. 

Rahmani et al. (2002) recently proposed a method for assigning dipole polarizabilities that 
takes into account the effects of target geometry, resulting in dipole polar izabilties that are location- 
dependent even when the target being approximated has homogeneous composition. Collinge & 
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Fig. 4. — As in Figure 1, but for m = 3 + 4z 

Draine (2004) have extended this treatment by adding finite-wavelength corrections to the "surface- 
corrected" polarizability prescription, and used DDA calculations for spheres and ellipsoids to show 
that accurate results for Q sca , and Q a b s can be obtained for modest numbers of dipoles even when 
the target material has a large refractive index \m\ > 5. The current calculations does not include 
such "surface corrections" , - the prescription (69) depends only on local properties - so it comes 
as no surprise that there are systematic errors even in the limit x — > 0. 

Collinge & Draine (2004) show that the error can be greatly reduced if the polarizabililites 
are assigned in a way that allows for both target geometry (i.e., surface corrections) and the finite 
wavelength corrections obtained from the lattice dispersion relation analysis. Using these "surface- 
corrected lattice dispersion relation polarizabilities" , Collinge & Draine (2004) show that for the 
case m = 3 + 4i, the error in (^abs can be reduced to of order 1% for x — 0.83 using as few as 
N = 5904 dipoles, whereas for x = 0.83, the error is seen from Figure 4 to be ~ 20% even with 
N = 17904 dipoles. Such surface corrections will obviously be of value for DDA calculations, but 
at the moment have been applied only to sperical, spheroidal, or ellipsoidal targets (Collinge & 
Draine 2004). 




1 1.5 2 2.5 

x = 27ra/A 



B. T. Draine :£004 03 01 



-14- 



Imlkd 

0.2 0.4 0.6 0.8 



1 

0.1 

0.01 
0.01 

0.005 



0.005 

0.01 
0.005 


0.005 
-0.01 

5 10 15 

X ~ 27Ta/A B.T. Dratne 2004.03.02 

Fig. 5. — As in Fi gure 2, but for N = 140408 dipoles. 

5. Summary 

The principal results of this paper are as follows: 

1. We have derived the exact mode equation (31) for propagation of electromagnetic waves 
through a rectangular lattice of polarizable points. 

2. If the electric polarizability tensor a for the lattice points is specified, then we can solve the 
mode equation to obtain the dispersion relation for waves propagating on the lattice. Alter- 
natively, if we wish the lattice to have the same dispersion relation as a continuum medium 
of complex refractive index m, then we can solve the mode equation for the appropriate 
polarizability tensor. 

3. In the case of a cubic lattice, our derived "lattice dispersion relation" polarizability tensor 
(69) is diagonal but anisotropic, and differs from the isotropic polarizabilities obtained by 
DG93. 

4. The new polarizability prescription has been tested in DDA scattering calculations for spheres. 
The resulting cross sections have accuracies comparable to those obtained with the DG93 
prescription. 
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A. Numerical Results for Ro(i), Ri, R2(i), and R^(i,j) 

The coefficients Rq{i) are obtained numerically from the sums (45). Note that eq. (45) for 
i?o involves the difference between two sums which separately diverge at infinity. To handle this 
divergence, one notes that in nature one is always interested in finite wavelengths and nonzero 
frequencies; when dealing with finite wavelengths, the lattice sums converge because of the oscil- 
lations introduced at large distances by retardation effects. To obtain the proper behavior when 
going to the zero frequency limit, we modify the summations to suppress the contribution from 
large distances. The manner in which this is achieved is not critical. We compute 



Ro(i) = lim 



E^-p(-«h 4 )-E||2-p(-«iqi 4 ) 



(Al) 



The exponential cutoff allows one to limit the summations to a|n| 4 < 20 and a|Q| 4 < 20. The 
summation (Al) is carried out for various small values of a and it is found, as anticipated, that 
for a < 10 -8 the numerical result for Rq is insensitive to the choice of a. Each of the sums in 
(Al) is over a lattice; |n| 2 defines the distance from the origin for the cubic lattice, and |Q| 2 is 
the distance from the origin for the reciprocal lattice. Each sum diverges at infinity, where we can 
replace the discrete summation by an integration; the reason why the divergences cancel is that 
the two lattices have the same density of lattice points, so that at large distances, both approach 
the same integration. 



Table 1: Ro(i) for selected values of d± : di : cfe 



dl 


d 2 


d 3 


Ro(l) 


R (2) 


Ro(3) 


1 


1 


1 











1 


1 


1.5 


0.20426 


0.20426 


-0.40851 


1 


1.5 


1.5 


0.52383 


-0.26192 


-0.26192 


1 


1 


2 


0.38545 


0.38545 


-0.77090 


1 


1.5 


2 


0.81199 


-0.21028 


-0.60172 


1 


2 


2 


1.19693 


-0.59846 


-0.59846 


1 


1 


3 


0.74498 


0.74498 


-1.48995 


1 


1.5 


3 


1.38481 


-0.14304 


-1.24176 


1 


2 


3 


1.96224 


-0.69714 


-1.26510 


1 


3 


3 


3.11030 


-1.55515 


-1.55515 



The quantities R\, R2, and R3 appearing in the dispersion relation in the limit of small but 
nonzero frequency v have been obtained by direct evaluation of R^{i^j) using eq. (49) followed 
by use of equations (50) and (51) to obtain R\ and i^W- Six elements of the symmetric matrix 
Rs(i,j) are given in Table 2. 
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Table 2: R 3 (i,j) for selected values of d\ : d 2 : d 3 



d\ 


d 2 


^3 


R 3 (l, 1) 


#3(2,2) 


#3(3,3) 


#3(1,2) 


#3(1,3) 


#3(2,3) 


1 


1 


1 


0.00000 


0.00000 


0.00000 


0.00000 


0.00000 


0.00000 


1 


1 


1.5 


0.37743 


0.37743 


-1.62922 


0.13566 


0.01609 


0.01609 


1 


1.5 


1.5 


0.80161 


-0.80815 


-0.80815 


0.16648 


0.16648 


-0.18041 


1 


1 


2 


0.55693 


0.55693 


-3.72878 


0.23735 


0.09360 


0.09360 


1 


1.5 


2 


1.02166 


-0.55501 


-2.30887 


0.27206 


0.25987 


-0.21806 


1 


2 


2 


1.26456 


-1.78732 


-1.78732 


0.37528 


0.37528 


-0.39535 


1 


1 


3 


0.81662 


0.81662 


-8.83412 


0.39793 


0.23223 


0.23223 


1 


1.5 


3 


1.34832 


-0.40088 


-6.06792 


0.43771 


0.42272 


-0.15845 


1 


2 


3 


1.62638 


-1.56624 


-4.80661 


0.55590 


0.55480 


-0.48211 


1 


3 


3 


2.04073 


-3.94766 


-3.94766 


0.76142 


0.76142 


-0.90991 



Finally, R\ and R 2 {i) are given in Table 3. 
Table 3: R\ and R 2 (i) for selected values of d\ : d 2 : d 3 



d 1 


d 2 


d 3 


#1 


#2(1) 


# 2 (2) 


#2(3) 


1 


1 


1 


0.00000 


0.00000 


0.00000 


0.00000 


1 


1 


1.5 


-0.53869 


0.52918 


0.52918 


-1.59705 


1 


1.5 


1.5 


-0.50962 


1.13457 


-0.82209 


-0.82209 


1 


1 


2 


-1.76582 


0.88788 


0.88788 


-3.54158 


1 


1.5 


2 


-1.21448 


1.55359 


-0.50100 


-2.26706 


1 


2 


2 


-1.59967 


2.01512 


-1.80739 


-1.80739 


1 


1 


3 


-5.47612 


1.44677 


1.44677 


-8.36967 


1 


1.5 


3 


-3.71651 


2.20875 


-0.12162 


-5.80365 


1 


2 


3 


-3.48931 


2.73708 


-1.49246 


-4.73393 


1 


3 


3 


-4.62875 


3.56356 


-4.09616 


-4.09616 
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